Clustering of inelastic soft spheres in homogeneous turbulence 
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In this paper we numerically investigate the influence of dissipation during particle collisions in an 
homogeneous turbulent velocity field by coupling a discrete element method to a Lattice-Boltzmann 
simulation with spectral forcing. We show that even at moderate particle volume fractions the 
influence of dissipative collisions is important. We also investigate the transition from a regime where 
the turbulent velocity field significantly influences the spatial distribution of particles to a regime 
where the distribution is mainly influenced by particle collisions. 
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I. INTRODUCTION 

Particles suspended in a fluid can show rather complex 
behavior and an accurate description of these effects is 
a longstanding challenge. In homogeneous isotropic tur- 
bulence preferential concentration [IH1] is a well known 
effect. There small, heavy particles tend to concentrate 
in regions where the vorticity of the fluid velocity field is 
low and the strain is high. This effect is most pronounced 
when the Stokes number St of the particles is St « 1. In 
the case of turbulence the Stokes number [5] depends on 
the particle and fluid densities p p and pf as well as the 
particle radius R p and the Kolmogorov length scale rj: 
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18 Pi 



Rr 



(1) 



Preferential concentration has been investigated inten- 
sively by experiments [TJ [2] and numerical simulations 
[BHS]. The occurrence of preferential concentration is 
closely related to the difference in inertia between the 
particles and the fluid, which is expressed by the Stokes 
number ( IJ . For small St particles follow the fluid stream 
lines rather closely and therefore, at least in an incom- 
pressible fluid, initially homogeneous distributed particles 
remain distributed that way. For large Stokes numbers 
the particle motion is only weakly influenced by the fluid 
which then leads to a diffusion-like motion of the particles. 
For intermediate Stokes numbers the local structure of 
the fluid velocity field becomes important |7 . In the case 
of dilute particle suspensions, when collisions between in- 
dividual particles can be ignored, heavy particles tend to 
concentrate in regions of high strain rate and low vorticity, 
and light particles tend to do the opposite The mech- 
anism behind this process is that the vortices present 
in a turbulent flow act as centrifuges that eject heavy 
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particles and attract light ones This effect is most 
pronounced for St « 1 pj [2j 14] . In most simulations one 
concentrates on an accurate description of the fluid and 
particle collisions are either ignored completely or treated 
as perfectly elastic. But in fact inelastic collisions are 
another mechanism that can lead to a clustering of par- 
ticles - an effect known as collisional cooling. As shown 
by Luding et al. jTQl [TT] for freely moving particles, this 
effect of (free) cooling is already important for moderate 
amounts of dissipation. 

In this paper we want to investigate the influence of 
dissipative collisions on the clustering of soft spheres in 
homogeneous turbulence. Additionally the dependence on 
the particle volume fraction shall be considered. We start 
by presenting the coupling of a discrete element model to 
a Lattice-Boltzmann simulation where a spectral forcing 
technique P21 [32] is used to generate turbulence. We 
then investigate the influence of the dissipation during 
collision on the clustering of particles for different Stokes 
numbers. Finally we vary the volume fraction by changing 
the number of particles in our system and study its effects. 



II. MODEL DESCRIPTION 

The discrete element model (DEM) [33] is a widely 
used and well established method for simulating granular 
materials 15 171. We use this model here to evolve a set 



of rip spherical particles with positions , masses m; and 
radii Ri according to Newtons' equation of motion 



(2) 



where is the total force acting on particle i. This force 
is given by the sum of a collision force and a drag force 
due to a fluid 



F . = pcoll + F arag_ ^ 

The collision force is given by a sum of two-particle colli- 
sions 
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and these two-particle collisions are modeled by a linear 
spring dash-pot model. In this model one first calculates 
the overlap between two particles i and j 



quadratically depend on the velocity difference. A widely 
used and well established empirical drag law for the tur- 
bulent case [TBH2"0] is given by 



5ij = (Ri + Rj) - |xj - Xj| 
and if Sij is positive, a repulsive dissipative force 



F C ° 11 = (fc n 



= ( kaSij + CnSij ) n y 



(5) 
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in the direction 



where the drag coefficient Cd is given by 
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acts on particle i. The model corresponds to a damped 
linear spring with stiffness k n and damping coefficient c n . 
The value 8^ is given by the relative velocities of the two 
particles in the direction of ny, i.e. 



(v< -vj)- n 
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where Vj is the velocity of particle i. Instead of the 
damping coefficient c n it may be more convenient to work 
with the coefficient of restitution e n , which is an easier to 
handle material parameter. This coefficient is a measure 
for how much energy is retained after a collision. To relate 
c n with e„ we additionally introduce the collision time t c 
as 



t c = — ( 7r — 2 arctan — 



(9) 



where 



with 



rriij 4m?- 



in. 



and r? c = (10) 
Irriij 



(11) 



rrii + rrij 

The coefficient of restitution is then given by 



(12) 



For the drag force F^ rag we use an empirical drag law. 
For a laminar flow the drag force [TH] should be propor- 
tional to the difference between particle velocity and 
fluid velocity u 



U - V, 



where ti is the particle response time 

18 Vf pf 



(13) 
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Here Vf is the kinematic viscosity of the fluid. For a 
turbulent flow on the other hand the drag force should 



This coefficient depends on the particle Reynolds number 
Re p , which is given by 



Re p = 



2Ri\u- vt 



(17) 



Such a drag law describes a simplified dynamics for the 
particles, where the added mass effect as well as the 
Basset-Boussinesq history force are neglected. Such an 
assumption is reasonable for heavy particles. There is 
also no pressure gradient term present, since the fluid 
field is incompressible. 

To calculate the fluid velocity u we use a Lattice- 
Boltzmann method with spectral forcing O |T3]. In 
recent years the Lattice-Boltzmann (LB) method has be- 
come a very successful method for solving many different 
problems in fluid dynamics [2"lTf2"3"] . The standard LB 
equation with an external force g (x, t) reads 

/ i (x + c i *t J t + <yt)-/ i (x,t) = 

- u 6t(fi (x, t) - fP (x, t)) + Gi (x, t) , (18) 

where /, (x, t) is the probability of finding a particle at 
time t and lattice site x, that is moving in the direction of 
the discrete lattice velocity c^. Here i = 1, . . . , Q indexes 
the discretization of velocity space. Equation (181 is a 
discretized version of the Boltzmann-equation where the 
left-hand side relates to the free-streaming and the right- 
hand side is an approximation of the collision operator 
plus a forcing term Gi. Here the well established Bhatna- 
gar-Gross-Krook relation [24] was used, which models a 
simple relaxation of the population fi toward a local equi- 
librium /? q with relaxation frequency lo. The equilibrium 
distributions are given by a Maxwell-Boltzmann distri- 
bution which is expanded up to second order in terms of 
Hermite polynomials resulting in 



fi = UiPf 1 
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where u>i are lattice weights and c s is the speed of sound. 
In this paper we chose a D3Q19 lattice depicted in Fig. [l] 
The lattice weights for this lattice are given as 
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FIG. 1. Definition of the lattice velocity vectors ci, 
for the D3Q19 lattice. 
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and the speed of sound is 
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The fluid density pf and the velocity u are related to the 
probability distributions fi through [25] 



PfU = ^/;C 4 + — g 



(22) 
(23) 



The choice for the forcing terms Gi is not trivial, since one 
has to construct Q values Gi from a three-dimensional 
vector g. The choice should guarantee, that the correct in- 
compressible Navier-Stokes equations are recovered when 
performing the Chapman-Enskog expansion. Here we are 
using the expression by Guo et al. [25] 



r< / -| 1 \ fCi - U (Ci • u) 
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(24) 



Finally we have to specify the external force g, which 
is the driving mechanism for the fluid. This force has to 
be chosen such that an homogeneous and isotropic tur- 
bulent velocity field is generated. This task is frequently 
encountered in direct numerical simulations (DNS) and 
other cases (see e.g. [12l [26] and references therein). In 
this work we use a method introduced by Alvelius |12j 
which has already been used with the LB method by ten 
Cate et al. [T3] . Here only a brief review of the technique 
is given. 

The basic idea for calculating g (x, t) is to generate 
in Fourier space a random, divergence free force field 
g (k, t) that is only active at small wave- vectors k, and 
then take the inverse Fourier transform to get the force in 
real space. This corresponds to the picture of turbulence 
where energy is injected into the system at small k, which 



is than transported to higher and higher frequencies until 
it is finally dissipated. To fulfill the condition 



kg = 



(25) 



we have to choose a force g which is always perpendicular 
to k. This can be achieved by choosing 

g (k, t) = A ran (k, t) ei (k) + P ran (k, t) e 2 (k) , (26) 

where A lan and P ran are two random amplitudes and the 
unit vectors ei and e 2 are chosen perpendicular to k and 
each other. Alvelius [H] choice for these vectors is 
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and 



e 2 = 




(27b) 



where k = |k|. The amplitudes A ran and B ran are chosen 
as 
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exp(i0i)sin(</>) (28a) 



exp (i0 2 ) cos ((f)) , (28b) 



where <j>, 0i, 2 G [0, 27r] are chosen at every time step as 
uniformly distributed random numbers. The spectrum 
function S (k) should only be active in an interval [k a , fcb] 
at small wave number and is given by a Gaussian 



S(k) 



'Aexp(-fcMT) ke[K,k h ] 



(29) 



where c determines the width of S (fc), kf the position of 
its maximum, and 
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The value P m specifies the mean power input by the 
spectral force. The total power input during one time 
step actually consists of two contributions 



■ftotal = \ (g{ti) 'g(*i)> 



(g(t i )-u(t j )) = Pi + P 2 . 



(31) 

The first term is due to the force-force correlation and the 
second one is due to a force-velocity correlation. Since u 
is basically a random field, P 2 cannot be controlled and 
may become rather large. Therefore it is desirable to have 
P 2 = at every time step and consequently Pi = p n . 
This can be achieved by demanding 



u ■ g 







(32) 
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for every active wave- vector. The condition (32) can be 



fulfilled if the angles 0i and 02 in Eq. ( 28 ) are not chosen 
independently anymore. The angle 9\ is then given by 



tan 0i = 



sin cf> Re £i + cos <fi (sin rp Im £2 + cos ip Re £2) 



sin <j) Im £1 + cos (j> (sin ip Re £2 ~ cos if) Im £ 2 ) ' 

(33) 

where £1 = u • ei, £2 = u ■ 62, and ^1 € [0, 27r] is a 
uniformly distributed random number. The angle 02 is 
then calculated as 2 = ip — Q\ ■ 

By specifying three input parameters we can determine 
all the properties of the spectral force g: 

• The length l g , which defines the large scales of the 
turbulent field. 

• The characteristic velocity u g , which fixes the 
timcscalc of the simulations. 

• The Kolmogorov scale rj which determines the small- 
est turbulence scale. 

To ensure that the simulated fluid field is incompressible, 
the Mach number has to be small. This gives the condition 
that the characteristic velocity has to be much smaller 
than the speed of sound, i.e. u g <C c s . From l g we can 
determine fcf = n^/lg, where is the linear size of the LB 
lattice, and further fc a = 1 and fcb = 2 fcf — 1. In Eq. (29) 
the constant c is chosen as c = 1. The power input Pj n 
is specified by P xn — u g /l g and the fluid viscosity is then 
given by 
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(34) 



From the viscosity the relaxation frequency in Eq. ( 18 1 is 
determined as 
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(35) 



At this point it is worth noting that there is no feedback 
from the particles to the fluid. According to Elghobashi 
[27] a four-way coupling would be needed in the range of 
Stokes numbers and volume fractions considered in this 
paper. Unfortunately incorporating a feedback force into 
the simulations is not trivial. Adding the negative drag 
force that acts on each of the particles to the spectral 
force G would be one way of coupling the fluid to the 
particles. The resulting forces would be singular and not 
necessarily located on a specific site of the LB lattice. 
Nash et al. [25] used a regularized Dirac delta function 
to incorporate singular forces into the LB method in the 
case of low Reynolds numbers. The spectral force used 
in this paper is random and supplies energy to the whole 
system (on average) homogeneously distributed in space. 
The origin of the spectral force is completely artificial and 
energy conservation is therefore only fulfilled on average 
in time. Additionally, since the fluid velocity field is a 
random vector, the drag force and hence the feedback 



force are random vectors as well. Therefore including 
the feedback from the particles to the fluid would add an 
additional random component to the spectral force which 
could lead to an enhancement or a reduction of the local 
turbulent intensity. We assume that these effects cancel 
each other on average an thus the lack of the feedback 
force should not influence our results much. 



III. SIMULATION AND RESULTS 

We used the method described in the last section with 
a lattice of linear size til = 128. The forcing length scale 
was l g = 32, the characteristic fluid velocity u g = 0.05 and 
the Kolmogorov length scale 77 = 0.5. These values are the 
same as used by ten Cate et al. 
results in a fluid viscosity V{ 



. Using Eq. ( 34j) this 
0.00625. The fluid density 
was set to pf — 1.0 and the initial fluid velocity was u = 0. 
After initialization we let the LB part of the simulation 
run until the mean kinetic energy of the fluid stayed 
constant over a certain time. Into this turbulent velocity 
field we put n p particles of radius R p = 77 on a regular 
grid. The initial velocities of the particles were chosen 
equal to the local fluid velocity. Since the particles can be 
located anywhere inside the system and are not confined 
to the lattice sites of the LB grid, the local fluid velocity 
was calculated by linear interpolation from the velocities 
of the eight surrounding lattice points. We simulated 



four different number of particles 



TIt, 



884736, 262144, 



110592, and 32768) corresponding to volume fractions of 
(f) vol = 22.1%, 6.5%, 2.8%, and 0.8%. Additionally the 
Stokes number of the particles was varied to be St = 0.32, 
0.56, 1.0, 1.78, 3.16, 5.62, 10.0, 31.62, and 100.0. Using 
Eq. ([!]) it is possible to calculate the corresponding particle 
densities. The spring stiffness was set to fc n = 25.0 and 
the coefficient of restitution was varied as e n = 1.0, 0.95, 
0.60, and 0.25. Additionally we performed simulations 
where particle collisions have been ignored and therefore 
particles could overlap and cross each other. To quantify 
the clustering of particles we follow the work of Fessler 
et al. [2 . The main idea of this analysis is to measure 
the local deviation of the particle number from a Poisson 
distribution. This is reasonable since if one divides a 
cubic volume where particles are uniformly distributed in 
space in smaller cubes, the number of particles in these 
smaller boxes are Poisson distributed. In more detail we 
take a snapshot of the system, divide it into m 3 boxes of 
linear size l rn = n^/m, (to = 4, . . . , n^) and then count 
the number of particles s in each box. For every box size 
we then measure the mean /z m , the standard deviation 
<r m , as well as the distribution LI m (s) of the number of 
particles. If the particles are uniformly distributed in 
space, II m (s) should be a Poisson distribution 



IR 
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X s e~ 



(36) 



Fig. [2] shows two examples of these distributions. In 
Fig. |2T(a) the measured distribution is close to the Pois- 



5 




FIG. 2. Histograms of the number of particles for different Stokes number and different box size. The left figure shows a case 
where the particle number distribution is close to a Poisson distribution which means, that particles are uniformly distributed 
in space. On the other hand the right figure shows a case where the distribution is much broader and differs strongly from a 
Poisson distribution. In this case preferential concentration occurs. 



son distribution and therefore the particles are almost 
uniformly distributed in space, i.e. no preferential con- 
centration is visible. In Fig. [2] (b) on the other hand the 
measured distribution is much broader than a Poisson 
distribution with the same mean value \i m . This means 
that there are boxes with too large or too small numbers 
of particles to be compatible with a Poisson distribution. 
Therefore particles are clustered in small regions of space, 
i.e. preferential concentration occurs. To further quantify 
the deviation from a Poisson distribution it is useful to 
define [5] the value 

D m = ° m - <Tpoisson . (37) 

Mm 

Depending on the value of D m one can distinguish three 
regimes: 

D m > 0: The distribution Tl m (s) is broader than 
Ilpoisson (s; Mm) and therefore clustering occurs. The 
larger D m is the more pronounced this effect is. 

D m = 0: The distribution H rn (s) and IIp i sson (s; Mm) are 
equally broad. Particles are uniformly distributed 
in space. 

D m < 0: The distribution II m (s) is narrower than a cor- 
responding Poisson distribution. This means that 
all boxes contain more or less the same number of 
particles. 

We first applied this analysis to the case where colli- 
sions between particles were ignored. The behavior of 
D m for particles with different Stokes numbers St at a 
volume fraction of 4> vo \ = 22.1% is shown in Fig. [3] (a). 
The deviation of the particle distribution from a Poisson 
distribution is clearly visible and as expected [5] preferen- 
tial concentration is strongest at a Stokes numbers around 
St pa 1. 

In a next step we included collisions in the simulations. 
We first set the coefficient of restitution to e n = 1.0, which 



corresponds to elastic collisions. The results of this case 
are shown in Fig. [3] (b). The collisions reduce the strength 
of the preferential concentration, but the strongest cluster- 
ing is still visible for a Stokes number St « 1. This effect 
is to be expected since collisions introduce a mechanism 
that tries to move particles apart and counteracts the 
"attraction" of particles in regions of low vorticity and 
high strain. Therefore the overall strength of preferential 
concentration is reduced. Additionally D m even becomes 
negative for small box sizes n^/m. In contrast to the 
collisonless case, particles repel each other and cannot 
come arbitrarily close together. For box sizes close to 
the particle diameter this means that fluctuations in the 
particle number are strongly reduced. Therefore D m can 
become negative for these box sizes. 

Now we reduced the coefficient of restitution. As shown 
in Refs. [10l [11] the influence of the dissipation during 
particle collisions should already influence the clustering 
at moderate values of e n . Therefore we chose e n = 0.95. 
The results of these simulations are shown in Fig. [3] (c). 
One can see that the dependence of the clustering on 
the Stokes number is less pronounced in this case, since 
the curves of D m are rather close to each other for dif- 
ferent Stokes numbers. Even so maximal preferential 
concentration at St sa 1 can still be observed. 

Fig. [3] (d) finally shows the results of simulations where 
the coefficient of restitution was set to e n = 0.25. Here the 
maximal preferential concentration at St w 1 cannot be 
observed anymore, and the clustering is clearly dominated 
by particle collisions. The effect of clustering increases 
with increasing Stokes number, since particles with higher 
St are less influenced by the fluid velocity field. 

To further investigate the two regimes of preferential 
concentration due to turbulence and clustering due to 
collisional cooling we again followed Ref. [2] and deter- 
mined the maximum of D m . In Fig. [4] we then plot the 
values of D max versus St for different coefficients of resti- 
tution and also the collisionless case. In the latter case 
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FIG. 3. Dependence of D m on the box size for different Stokes numbers St with and without collisions at a volume fraction 
4> vo i = 22.1%. Particles are randomly distributed if D m = 0. If D m > preferential concentration occurs and if D m < 0, on 
average, the same number of particles is located in boxes of corresponding size, (a) Simulations without collisions. Preferential 
concentration is maximal for St w 1. (b) Simulations with collisions and coefficient of restitution e n = 1.0, which corresponds to 
clastic collisions. Preferential concentration is reduced but the maximum around St w 1 is still visible, (c) Simulations with 
e n = 0.95. Clustering of particles is less dependent on the Stokes number, but the maximum of preferential concentration is still 
visible, (d) Simulations for e n = 0.25. Dissipative collisions dominate the clustering of particles. Since particles with higher 
Stokes number are less influenced by the turbulent fluid, clustering dominates. 



preferential concentration is most pronounced. This is 
clear since in this case there is no effect which tries to 
move particles apart. This means particle can overlap 
and therefore many particles can concentrate in regions 
of the fluid velocity field where vorticity is low and strain 
is high. When particle collisions are included in the simu- 
lation, large overlaps between particles are not allowed 
anymore. Therefore preferential concentration is reduced. 
This effect can also be seen in Fig. [4] for e n = 1.0. The 
maximum at St w 1 is still clearly visible but the overall 
strength of preferential concentration is reduced. When 
the coefficient of restitution is reduced, another effect for 
the clustering of particles is introduced. Therefore D max 
should increase when e n is reduced. This effect can be 
seen in Fig.|4]for e n = 0.95. The increase of D max is larger 
for higher Stokes numbers. This can be explained by the 
fact, that particles with larger St are less influenced by 
the fluid velocity field. Further decreasing e n brings us 



into a regime where the dissipative particle collisions are 
the dominant mechanism for clustering of particles. Fig. [4] 
shows the behavior of -D max for e n = 0.6 and e n = 0.25. 

Fig. [5] shows the same data as Fig. [4] but this time ex- 
hibiting the dependence of -D max on e n for different Stokes 
numbers. The plot shows a crossover at a coefficient of 
restitution around e n ;$ 0.95. Above this threshold the 
clustering is influenced by the turbulent velocity field 
and preferential concentration at St ~ 1 can be observed. 
Below the threshold the collisions between particles are 
the dominant mechanism for clustering of particles. 

We further investigated the influence of the particle 
volume fraction <f) vo \ on the clustering of particles. By 
reducing the number of particles n p we simulated systems 
with volume fractions 4> vo \ — 6.4%, 2.8%, and 0.8%. Again 
we varied the Stokes number of the particles and measured 
Dm for different box sizes — . We then determined the 

m m 

maximum of D m and plotted the results in Fig. m\ The 
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collisions is still visible. At (j> vo \ = 6.5% with e n = 0.25 
the clustering is still dominated by the particle collisions. 
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FIG. 4. Plot of the maximum of D m for different Stokes 
numbers St at a volume fraction of </>voi = 22.1%. The simula- 
tion without collisions clearly shows the strongest preferential 
concentration at St ~ 1. This effect is reduced by collisions. 
For a coefficient of restitution of e n = 0.95 a maximum at 
St ~ 1 is still visible. For smaller e n this maximum disappears 
and the clustering becomes larger for larger Stokes number. 
The reason for this is that the fluid velocity field less and 
less influences the motion of the particles and the collisions 
become dominant. 
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FIG. 5. The dependence of D max on e n for different Stokes 
numbers. (Same data as in Fig. |4| At e n ~ 0.95 a crossover 
is visible. For e n above the crossover the influence of the 
turbulent velocity field is still visible and preferential concen- 
tration occurs; below, the main mechanism for clustering are 
the dissipative particle collisions. 



first thing to notice in these plots is the influence of the 
particle collisions. As expected the difference between the 
result of the simulations without collisions and the case 
of elastic collisions (e n = 1.0) becomes smaller with lower 
volume fractions. The difference between these two cases 
is less distinct for particles with larger Stokes numbers 
and more pronounced for St ss 1. For <p vo \ = 0.8%, the 
difference between these two cases even disappears for 
Stokes number above 10. 

The next point to mention is that even at moderate 
volume fractions the influence of the dissipative particle 



At the lower volume fraction 



= 2.8% the influence of 



the turbulent velocity field can be seen by a maximum of 
D max at St w 1, but the effect of the dissipative particle 
collisions is recognizable at larger Stokes numbers. Only at 
rather small volume fraction of (f> vo \ = 0.8% the influence 
of the particle collisions becomes negligible. 

As already explained in Section [XT] the model does not 
implement any feedback from the particles onto the fluid. 
The influence of the particles on the turbulent velocity 
field is known in the literature as turbulence modula- 
tion [3]. Depending on several parameters like particle 
radius, Stokes number, volume fraction, Reynolds number, 
etc., this modification may be either an enhancement or 
a reduction of the turbulence intensity, which in turn can 
influence the local structure of preferential concentration. 
Despite several experimental, numerical and theoretical 
investigations a conclusive understanding of the involved 
phenomena is still not found 4 . It is reasonable to as- 
sume that the modification is less important for lower 
volume fractions. Apart from that, it is rather difficult 
to reliably predict the changes due to turbulence modula- 
tion. Recent high-resolution simulations [22] and highly 
resolved particle image velocimetry experiments |30j show 
that for particles with diameter d gl r\ dissipation around 
the particles is strongly enhanced which then leads to a 
reduction in turbulence intensity in this regions. Includ- 
ing a feedback force in our simulations may therefore lead 
to a reduction of preferential concentration. 



IV. CONCLUSION 

In this paper we used a DEM together with a LB 
method to simulate the motion of inelastically colliding 
soft spheres in a homogeneous turbulent flow field. We 
investigated the influence of dissipative particle collisions 
on the clustering of particles and found that already at 
low densities collisions can be an important factor. For 
volume fractions around 22.1% collisions become domi- 
nant below a coefficient of restitution of e n ~ 0.95. Above 
this threshold we observed preferential concentration with 
a maximum for particles with a Stokes number around 
St w 1. For lower volume fractions the influence of par- 
ticle collisions becomes less important and the crossover 
between preferential concentration and collisional cooling 
is shifted to smaller values of e n . 
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FIG. 6. Dependence of -Dmax on the Stokes number St for different coefficients of restitution e n and the collisionless case at 
different volume fractions V ol- Even at moderate volume fraction the influence of the dissipative particle collisions is clearly 
visible. Only for densities below 1% the effect becomes smaller. Also visible in this figure is the difference between the results of 
the simulations with elastic collisions (e n = 1.0) and the collisionless case. 
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